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Abstract 

An important issue in theoretical epidemiology is the epidemic threshold 
phenomenon, which specify the conditions for an epidemic to grow or die out. 
In standard (mean-field-like) compartmental models the concept of the basic 
reproductive number, Rq, has been systematically employed as a predictor for 
epidemic spread and as an analytical tool to study the threshold conditions. 
Despite the importance of this quantity, there are no general formulation of 
Rq when one considers the spread of a disease in a generic finite population, 
involving, for instance, arbitrary topology of inter-individual interactions and 
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heterogeneous mixing of susceptible and immune individuals. The goal of this 
work is to study this concept in a generalized stochastic system described in 
terms of global and local variables. In particular, the dependence of Rq on 
the space of parameters that define the model is investigated; it is found that 
near of the "classical" epidemic threshold transition the uncertainty about 
the strength of the epidemic process still is significantly large. The forecast- 
ing attributes of Rq for a discrete finite system is discussed and generalized; 
in particular, it is shown that, for a discrete finite system, the pretentious 
predictive power of Rq is significantly reduced. 
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I. INTRODUCTION 



It is nowadays recognized that the phenomenon of health- disease in human communities 
only may be understood by considering complex and dynamic inter-relations among several 
factors operating simultaneously in multiple spatiotemporal and organizational scales. In 
fact, the healthy and sick individual suffers uninterruptedly the effects of the microbiological 
evolution, the antropogenic environmental and ecosystem stress and many others misdeeds 
resulting from socioeconomic inequalities. Therefore, it is not surprising to find out the 
proliferation of a myriad of methodological tools employed during the development of the 
epidemiological research. 

Among this methodological mosaic the mathematical and computer (or simulation) mod- 
eling of communicable and infectious disease comes as a hypothetical-deductive approach 
whose scope consists primarily in understanding and manipulating, a priori and to pre- 
dictive purposes, the underlying mechanisms behind the origin and diffusion of epidemic 
events. As a matter of fact, the attempt of understanding in what conditions pathogenic 
agents (once invaded a host population) could establish themselves as an infection (the 
transmission of pathogens from one host to another) resulted in the development of one 
of the most important and thoroughly discussed concepts in infectious disease modeling 
as early as in the beginning of the last century, namely the epidemic threshold. Thus, in 
writings of R. Ross (1909) [1] the so-called mosquito theorem was the first recognition of a 
quantitative threshold deducing that it was not necessary to eliminate mosquitoes totally 
in order to eradicate malaria. Two decades later would testify the publication of the classic 
Kermack-McKendrick's (1927) [2] paper that definitely consolidated the threshold concept 
in epidemiologic literature. In this deterministic SIR model [S stands for susceptibles, / 
for infected, and it! for removed) an epidemic process is considered to evolve only when the 
density of susceptible individuals is greater than a threshold value Sc- Bartlett (1957) [3], 
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based in a large amount of collected data of disease incidence in industrialized countries 
introduced thirty years later another expression linking microbial invasion and threshold 
parameters: the critical community size, that could explain the fade-out patterns of measles 
epidemics. 

However, the inherent individual heterogeneity and probabilistic local nature of interindi- 
vidual relationships has been traditionally neglected in state- variable models like this; in fact, 
in this population level approach all behavioral and individual variability are diluted into 
the intercompartmental rates and densities or number of mean individuals — as S, I or R 
compartments — described in terms of partial or ordinary differential equations. Neverthe- 
less, it was subsequently possible to express the epidemic threshold in a way perhaps much 
more intuitive when the focus changed to consider the infected host or the parasite itself, 
instead of looking at the density or number of susceptible. In this perspective the thresh- 
old condition that determines whether an infectious disease will spread in a susceptible 
population has been described through the so-called basic reproduction number or also de- 
nominated as basic reproductive rate, commonly denoted by TZq [4] . For microparasites such 
as viruses or bacteria it may be biologically understood as the average number of secondary 
cases produced or caused by one infected individual during its entire infectious period in 
a completely susceptible population. Thus, the intrinsically individual based perspective 
of this threshold concept should not be underestimated since the reproduction number can 
link the inside-host evolutionary or pathogenic dynamic (microscale) and transmission pro- 
cess at population level (macroscale) . From a purely deterministic point of view it appears 
intuitively evident that if TZq > 1 the pathogen can undoubtedly establish itself in a host 
population and, at least, an endemic regime will settle down. But this is a short-sighted 
prediction since, specially to directly transmitted disease in finite populations, the mecha- 
nisms that ensure the maintenance of the parasite within a community depends critically 
on the way as the individuals interact one another, sometimes unforeseeable. 

In this work we analyze limitations of the predictive power of the Rq parameter (as classi- 
cally formulated) for the spread of a disease: Alternatively to population level approach and 
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state- variable models, stochastic inter-individual interactions are also used and its implica- 
tions on the predictive attributes of the basic reproduction number TZq are studied through a 
simplified model: a lattice based model including infectious period in that individual inter- 
actions are straightforwardly described in terms of global ( F ) and local (A) variables, which 
in turn can be tuned out to simulate respectively the populational mobility and geographical 
neighborhood contacts. 

The remainder of this paper is organized as follows. In the next section it is presented a 
general formalism to the evolution of a population invaded by an infection. The formalism 
is then applied in section 3, where concepts involving TZo and the threshold phenomenon 
are discussed in order to define an invasion criterion for the infection and evaluation of TZq. 
The results are discussed in section 4- Although this work will be mainly concerned on TZq 
as a function of the model's parameters, the formalism presented in what follows can be 
applied to study a variety of epidemic scenarios. 

II. THE MODEL SYSTEM 

Consider a discrete dynamical system (discrete space and discrete time) where a popu- 
lation of N individuals is distributed on the sites of a toroidal lattice M = {rriij} — ^with 
i and j varying from 1 to L {N — L x L). Each individual site rriij is assigned to re- 
ceive three personal specific attributes: (1) Si spatial address or lattice position (2) 
a set of three possible status, namely, s, i and r, specifying a clinic disease stage of each 
particular individual, which represent, respectively, the conditions of susceptible (subject to 
be infected by a contagious agent), infectious (effectively transmitter of contagious agents) 
and removed (recovered or immune); and finely (3) an infectious period r, specifying how 
many units of time an infected individual can propagate the contagious agent. Note that 
'^s + ^i + '^r = N, with N constant. 

Such a system is suitable mainly for describing a single epidemic in a closed system 
(no birth or migration). The choice of such reduced model, however, is not far-fetched 
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because, as already mentioned above, the main interest here involves only very short period 
of time, so that the dynamics of host births, migration, etc., are largely irrelevant. [5]. The 
dynamic evolution of the population is described, step-by-step, by a set of a priori stated 
interaction rules, and assumes that each new configurational state of the system (described 
here by the geographical address (i, j) of each individual and by the instantaneous number 
of susceptibles S{t), infectives /(t), and removed individuals R{t) ) depends only on its 
previous state. Hence, for the present purpose the spread of the disease in the population 
is considered as being governed by the following rules: 

1. Any susceptible individual may become infected with a probability ps- An infected 
susceptible becomes infective after an average latency time r; (assumed here as = 0, 
without lost of generality) . 

2. Infectives are removed deterministically from the system (becoming immune) after 
an infectious period r, that for simplicity is considered as constant for all infected 
individuals. 

3. Once in the removed class the individual participate only passively in the spreading 
of the infection (eventual topological blocking) by a period of immunity greater than 
the complete epidemic process. 

During one time step, the three preceding rules are applied synchronously to all sites in 
the lattice; the present model, therefore, can be viewed as a simple two-dimensional cellular 
automaton. Actually, it is an adaptation of automata network to standard SIR models for 
studying the spread of infectious diseases. 

In this work, the probability ps, which is intended to be probability per unit of time, is 
taken as the superposition of the local and global influences, in order to unify the individual- 
based (contacts among nearest neighbors) and the standard mean-field (homogeneously in- 
teracting population) approaches. Therefore, one assumes that disease transmission occurs 
with a total infection probability ps written as 
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PS = VpG + ApL, 



(1) 



where the pre-factors F and A are weight parameters tuning the short (cluster formation) 
and long-range (mean-field type) interactions; it is also required that F -|- A = 1 in order to 
satisfy the probabilistic requirement ^ ^ 1. 

The global influence pa amounts to the probability of a susceptible to become infective 
due to the ubiquity of I{t) infected individuals (mean- field). So one can expect that in the 
limit of large N [N — > oo), in each time step, any susceptible may become infected with 
probability 



where < p < 1 is one of the model parameters: it limits the maximum value of pc and 
is related to the intrinsic mobility of the population; the sum sweeps all lattice sites {k, I}, 
and 5i^cr{k,i) is the Kronecker delta function which assumes the value "one" when the state 
a of the site (A;, /) corresponds to the infectious state i, and "zero" otherwise (cr(A;, /) can 
be s,i or r). Actually, the sum in the Equation 2 just counts the instantaneous number of 
infectious individuals I{t) in the population. 

On the other hand, the local term pi = Phii-ij) is the probability of a susceptible 
individual (located at the site (i, j)) contracting infection due to n infectives first and second 
neighbors (0<n<8isa integer number corresponding to all possible combinations of 
+ ^) J + 0) with = 0, 1, —1). Therefore, let A G [0, 1] be the probability of a particular 
susceptible when just one of its neighbors is infective. Hence, (1 — A)" will be the probability 
for not contracting the disease when exposed to n infectives. Therefore, the chance of he 
(or she) contracting the disease in a unit of time is [6] 



Thus, when A = 1 the infection spreads deterministically, with 8 nearest neighbors to any 
infective being infected (the choice for equipotcnt first and second neighbors was adopted 
because the use of only the four nearest neighbors is unduly restricting). 




(2) 



{k,i} 



p^^l-{l-\f. 



(3) 
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The expression for pc is a convenient and simple way for describing the populational 
mobihty. It is based on the mass action law, borrowed from the chemistry, and gets new 
meaning here under the perspective of pairwise spatially disordered interactions through the 
population elements. In this sense, it is a result of the small-word effect, and so became a 
particular version of the small-word lattice of Watts and Strogatz [7] 

This simple approach allows to study in great detail the dynamical behavior of the model 
in the full space of control parameters A and p, and the local and global balance pre-factors A 
and r. Therefore, the system is governed by ps (Eq.l) and r, and its temporal evolution is 
determined by updating the lattice synchronously at each time step through the application 
of the three rules above. 

III. 7^o AND THE THRESHOLD PHENOMENON 

The probability pi as in the Eq.(3) [8], [6], and in a number of alternative forms [9], 
has been employed in the analogy between percolation and epidemic. Since that the critical 
value Pci in which random clusters grow to infinite size, is know (analytic or numerically) for 
any lattices, Pc may be used as a powerful general criterion for "epidemic spread" [10], [9]. 
However, due to the traditional importance of the concept of Rq in the epidemic scenario, 
this threshold is generalized for finite discrete systems, as described above, in order to show 
its relevance for an intrinsically individual based perspective of the problem. 

The overall structure of the model presented here shows the interplay of two types 
of transmission mechanisms by assuming that each infectious individual interacts strongly 
(physically) with their few susceptible neighbors, and uniformly and weakly with each par- 
ticular susceptible in the population of susceptibles. Thus, the local mode of transmission 
Pl incorporates the individual-based component from the perspective of the susceptible in- 
dividuals, the actual (physical) contacts that each susceptible experiences, and the global 
probability pg, due to intrinsic populational mobility, which may be viewed as resulting of a 
mean-field (discrete) approach, in the sense that the disease transmission to each susceptible 
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individual also depends on the instantaneous total number of infectious individuals in the 
population. 

To better appreciate the consequences of this formulation, it is firstly run simulations 
for the extreme values of the tuning pre-factors through the procedure described above. 
These cases allow to recover the two modes of transmission in its pure form corresponding 
to the {i) homogeneous mixing approximation (mean field), when F = 1, and to the {ii) 
percolation process (the transmission occurs by localized individual contact), when A = 1. 
Furthermore, it is considered the damage AI on a susceptible population due to just one 
infected individual (/(O) = 1) landing in a totally susceptible population {S{0) = N — 1 and 
-R(O) = 0) during the infectious period r; to calculate A7 it is considered only the number 
of new infected individuals in the population after r time steps, ignoring infections from the 
victims of these first infected individual (operationally, it is enough to consider the latency 
Ti > T, that is, a latent period of infection greater then the infectious period). The Figures 
la and lb show, respectively, the simulation results for A7 as function of p for F = 1, and 
the behavior of AI as a function of the contact probability A for A = 1 (that is, F = 0); the 
system size considered in most of the simulations presented here was L — 100 (population 
size N — L X L — 10^), although some extra different sizes L/2 and 2L) were also used in 
order to verify finite size effects. F = 1 and A = 1 are the two limiting cases usually taken as 
references in studying the effect on the system when both mechanisms are superposed; the 
amount A7 is obtained after r = 10 time steps (covering exactly the infectious period) and 
was estimated as an average over 31 independent simulations (what is equivalent to verify 
the establishment of infection on 31 distinct populations with the same pattern of contacts 
among the individuals). 

The linear pattern observed for AJ vs p means that the present stochastic approach 
reproduces qualitatively the classical basic reproductive number TZq if one identifies AI as 
the average number of secondary cases that an infectious individual causes. Indeed, the 
linear relation AI = [-^ S{0)] r fits pretty well the data shown in Figure la, and so one 
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may consider that infectives make contacts at a mean rate [-^ 'S'(O)] throughout an infectious 
period of length r (note that for large enough populations S{0) p). On the other hand, 
when A = 1, the amount AI represents TZq for the case where individuals interact only with 
their spatial nearest neighbors, and so its values saturates at A7 = 8 for A > 0.3. For 
each particular run, significant fluctuations on AI are observed (mainly for smaller A^) but 
averaging over 31 runs is enough to smooth considerably the curves, as shown in Figure 1. 

Before to proceed through the application of the present formulation, some comments re- 
garding the deflnition of the TZq are in order. The basic reproduction number has been widely 
used as a predictor parameter conceived to indicate the epidemic potential of a pathogen 
once it has introduced in a totally susceptible population. In fact, to deterministic and 
continuous (in space and in time) population-system models the future fate of an infectious 
agent has been expressed through the threshold condition. Accordingly, when TZq > 1, in- 
fections can invade a totally susceptible population and persist; if TZq < 1, the disease then 
dies out and can not establish itself. To the special condition TZq — 1, there is an endemic 
regime in that the unique initial infectious case reproduces subsequently just one infectious 
secondary case and son on. 

This assumption in modelling of the establishment of an infection (which is possibly 
wrong) [11] will be partially preserved here to have the classical treatment as a reference 
but, indeed, to capture more realistic or probable practical situations is of interesse that the 
"flrst analytical look" at a population be considered when the epidemic process is already 
in course. For instance, at the initial time to one may consider the arbitrary situation in 
that I{t — to) » 1 at the same time that the number of removed individuals is also large, 
and then ask the question: What is the value for the reproduction number in this case? To 
answer this question one may generahze the concept of TZq as the normalized average number 
TZ{to;T) of secondary cases (reproductive ratio) about the time due to I{to) infectious 
present in the population at t — to, through the following expression 





(4) 
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where the brackets means an average on a set of n independent runs in the time interval 
[to, to + ''"]; ^-iid the sum over {/c,/} sweeps all sites occupied by individuals in the status 
s (susceptibles) . Note that all the instantaneous extensive and intensive conditions of the 
population, at any arbitrary time to, are all taken into account, as for example, the sites in 
the removed status randomly scattered through the population (acting as epidemic shield 
protectors), and the set remaining infectious time r(A;, /; to) for each individual in the status 
i located at the site {k,l). These conditions certainly affect the epidemic process and the 
progression of the epidemic process depends in some how on the reproduction number's value, 
(that is, if TZ{to;T) > 1 or < 1). But, as already mentioned above, the initial condition 
I{t — to) — 1 will be deliberately used in the present work in order to maintain the original 
intention of comparing the traditional deterministic definition of the basic reproductive ratio 
TZq with the present stochastic approach. 

In order to infer how the intrinsic stochastic nature of the epidemic process affects the 
predictive attributes of TZq, the concept of epidemic probability Pe is introduced. Numerically 
it is estimated directly from the simulation experiments based on the algorithm of the 
previous section. Indeed, it is just given by the ratio Pe — rie/n, where rig is the number of 
runs in that at least one susceptibles was infected during the infectious period, and n is the 
total number of runs or experimental populations. The probability Pe may be expressed 
as function of the mean reproduction number TZq, which also is determined from the same 
simulation experiments by using the Equation 4 above. In the Figure 2 it is shown the 
resulting Pe as a function of TZq with F = A = 0.5 and p and A varying in the interval 
(0 — 0.2]. The large number of scattered points in the graph, mainly at larger TZq, is an 
intrinsic aspect of this graph due to the fact that in the parameter space (p. A) there are 
different combinations of p and A resulting in approximately the same values for TZq, as it is 
illustrated in Figure 3. Therefore only the stochastic component of such scattering of points 
may be reduced by increasing the number of runs used in the averaging procedure. 

The amount Pe tends to saturates at Pe — 1 when the value of TZq is sufficient large 
(T^-o > 3), so that the epidemic spread in the population almost always is observed. Fur- 
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thermore, the results showed in the Figure 2 means that only for large enough TZo (actually 
TZq > 3) one can be sure about an epidemic development in the population, while that, 
even for TZq < 1 there is still a possibility to have an epidemic spread. Therefore, from the 
epidemic control perspective, reducing the effective reproductive number to a level below 
one, upon vaccination, for instance, could be a potential problem of strategy since that for 
TZq < 1 in about 60% of events this strategy will fail, that is, an epidemic process should 
be estabhshed with chance of 60% for 7?.o = 1, under the conditions of the present model. 
More pointedly, despite the claim of the threshold criterium, it is improbable to recognize 
(using only standard census data) the imminence of any epidemic disaster if the system is 
near to the threshold region. [12] The more accurate (although frustrating) criterium is to 
realize that, irrespective the value of TZq that the level of vaccination forces, there is always 
a chance (even thought small) of the disease re-invading the population. 

The same system size N — 10"^ was employed in order to get all the results discussed 
above. However, in order to verify eventual effect of the system size on the results, two 
extra systems were considered, namely a smaller N = x 10"^ and a bigger N = 4 x 10^ 
system, but no significant difference was found. Clearly fluctuations are smaller for larger 
systems mainly because the chance of nucleation of closer cluster due to the global term Fpa 
decreases with the system size A^, reducing then the chance of the magnification effect of the 
local term Api on eventual clusters located nearly enough each other. The Figure 4a shows 
for Pe vs TZq (in the interval < T^-o < 2) for two different system sizes; note that the size 
effect is pronounced only on the second moment (dispersion of the data) of the distribution 
of Pe for each TZo- More precisely, the Figure 4b shows the normalized standard deviation 
(relative error) a^^ as function of TZq for the larger 4 x 10^ system. A decreasing I/TZq— like 
behavior for the relative error is a consequence of the averaging of integer quantities, that 
is: TZq — (0xno + lxni-|-2xn2 + 3x ^3+ ■■■)/■>]; where rj — no + ni+ 713+ • • •, and is 
the number of experiments in which exactly k susceptibles were infected. 

Finally, the numerical equivalence between TZ'^ estimated by an analytical approximation 
and TZo calculated by simulation is verified. For this purpose TZ'q is considered in the limit 
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of large populations {N — > oo) by taking the mean number of susceptible infected by just 
one infective during its infectious period r, through the following direct expression 



The Figure 5 shows the parametric graph of TZq vs TZq where they are calculated, respec- 
tively, by Equation (5) above and by simulation using the proposed probabilistic approach 
represented in Equation (4), with I {to = 0) = 1. Strong correlation between the two ways 
for estimating the basic reproduction number is kept only for values of p and A not too large 
(7^0 ^ 2) because during the time r, the local term that composes TZq (Eq.4) may change 
from zero up to eight, while this limit is not present in the Equation (5). However, that is 
enough in order to validate the conclusions about the predictive attributes of Pe — Pe{T^o) 
because TZ'q and TZq are numerically equivalent: the result given by Equation 5, although 
intuitive, follows from a stochastic representation of the classical SIR model [13]. 



In this paper a stochastic version of the original SIR model (involving only single epi- 
demics) was introduced with the main purpose of to characterize and re-interpret the con- 
ditions for the establishment of an epidemic in a population through the concept of basic 
reproduction number (TZq). A peculiar characteristic of the present approach is the assump- 
tion that the probability of a susceptible individual become infective is a superposition of the 
local and global influences. Using as initial configuration just one single infected individual 
in a fully susceptible population, condition frequently used to define TZq, it was demonstrated 
that the discrete character of a finite population reduces the pretentious predictability of the 
threshold criteria, and so it is, indeed, an incomplete predictive tool since that, irrespective 
to the value of TZq, an epidemic has a finite probability to establish itself, due the inherent 
stochastic nature of any finite epidemic system. 




(5) 



IV. FINAL COMMENTS 
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Indeed, more consistent derivation of Rq has been tried, even though using the same 
classical deterministic approach, due to the too widely estimate obtained to Rq, which in 
recent applications for the smallpox have varied from i?o = 1-5 to > 20 [14]. Rather than 
just a caricature of the original formulation of TZq, the approach presented in this paper 
may be viewed as a simpler and generic alternative for investigating the spread of diseases 
in a population, which may greatly facilitate the analysis of a number of distinct epidemic 
scenarios. Particularly, a system with increasing topological complexity can be easily tackled. 
For example, one may consider the practical situation in that, at an arbitrary initial time 
to, the population has already many infectious individuals (that is, /(to) >> 1), and also 
many immunes scattered through the population (working as epidemic shield) and then try 
to answer the question: What is the value for the reproduction number in this case? 

Finally, as a major challenge that this "microscopic" approach can handle, one may 
think on the possibility of incorporating in the traditional definition of TZq the underlying 
evolutionary dynamics of the pathogenic agent. This view is in contrast with the standard 
epidemiological models, which tend to use a constant absolute parasite fitness TZq. However, 
more detailed considerations on the investigation of this avenue of research is left for a future 
contribution. 
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FIGURE CAPTIONS 



FIGURE 1: Average damage A7 due to just one infectious individual on the susceptible 
population S = N — 1, ior two extreme cases . [a]- F = 1: the amount A/ changes hnearly 
with the intrinsic mobility p, as can be expect from Equation (2). [b]- A = 1: the amount 
A7 increases rapidly with the infection probability A due to local (physical) contact, and 
saturates at A7 = 8 for A > 0.3, as one can infer from Equation (3). 

FIGURE 2: The epidemic probability as a function of average reproduction number TZq. 
The tuning pre-f actor are fixed at F = A = 0.5, and the parameter p and A are choosing 
from the interval [0, 0.2]. For T^-o ~ 1 epidemics are observed in about 60% of the events (in 
a population of size N — lO'^). 

FIGURE 3 Reproduction number TZq as function of the model parameters (p, A) obtained 
by averaging over 100 independent realizations. Each strip, identified by a different gray 
tone, corresponds to a range of value for TZq according to: white, < T^-o < 1; hght gray, 
1 < 7?.o < 2; and so on. At the limit of very large populations {N — > oo) the slope a (dotted 
lines), which roughly delimitates each region, can be obtained using Equation (5) -see text; 
giving a = ^. Therefore, in [a] F = A = 0.5, so a = -0.125; and in [b] F = 0.9 and A = 0.1, 
giving a — —1.125, whose values are closely reproduced by the results. 

FIGURE 4: [a] The epidemic probabihty Pe vs TZq for two systems: N = A"'^ x 10^ (open 
circles) and N = 4x 10^ (dark circles) smaller fluctuations for the larger system is the most 
significative difference. [&]- The relative error decreases as TZq increases; for 7?.o ~ 1 the 
absolute error is of the same magnitude of TZq as a consequence of the averaging on "zeros" 
and "ones", mainly. 

FIGURE 5: Numerical equivalence between TZq estimated by an analytical approximation 
and TZq calculated by simulation. 
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